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Abstract 

Clean superconductors with weakly coupled conducting planes have 
been suggested as promising candidates for observing the Fulde-Ferrell- 
Larkin-Ovchinnikov (FFLO) state. We consider here a layered super- 
conductor in a magnetic field of arbitrary orientation with respect 
to the conducting plane. In this case there is competition of Pauli 
spin-pair-breaking effects, favoring the FFLO state, and orbital-pair- 
breaking effects, favouring the Abrikosov vortex state. In previous 
work, phase transitions to phases with pairing in Landau levels with 
quantum numbers n > have been predicted. Here, we calculate 
the actual structure of the stable states below H C 2 by minimizing the 
free energy. We find new order parameter structures differing from 
both the traditional Abrikosov and FFLO solutions. These include 
two-dimensional periodic structures with several zeros of the order 
parameter, as well as quasi-one-dimensional structures consisting of 
vortex chains separated by FFLO domains. We discuss the limit of 
high n, where some interesting but yet unsolved questions appear. 

PACS numbers: 74.25.Ha, 74.80.-g, 74.80.Dm 
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1 Introduction 



Suppression of superconductivity by a magnetic field is a consequence of its 
interaction with either the magnetic moment due to the orbital angular mo- 
mentum J]] or due to the spin of the electrons. The recent discovery of 
several classes of layered compounds with nearly decoupled superconducting 
planes has renewed the interest in the interplay of spin and orbital pair break- 
ing mechanisms. These include high-T c cuprate superconductors, organic su- 
perconductors with very large upper critical fields, such as (TMTSF) 2 PF 6 0, 
and hybrid ruthenate-cuprate compounds like RuSr 2 GdCu20s @, ||. The or- 
bital pair breaking effect, described by a vector potential A, is usually much 
larger than the spin effect and dominates most of the phenomena observed 
in magnetic fields, e.g. the structure of vortices in type II superconductors. 
Nevertheless, the spin effect has attracted considerable interest 0, |J; it has 
soon be realized that its relative importance can be greatly enhanced by re- 
ducing the spatial extension of samples in a direction perpendicular to the 
external field. For the case of an almost two-dimensional superconductor 
with applied field parallel to the conducting planes, the orbital upper critical 
field is extremely high and the spin pair-breaking mechanism becomes the 
dominant one. 

Then, if the orbital pair breaking mechanism can be completely neglected, 
the homogeneous superconducting state at T = becomes energetically 
favourable 0, |6| at an upper critical field given by fxH c = A /v^2, where 
A is the BCS gap at T = (Chandrasekhar-Clogston limit). Later, a 
spatially inhomogeneous state (FFLO state) with higher critical field was 
theoretically predicted by Fulde and Ferrell 0, and Larkin and Ovchin- 
nikov 0]. In the FFLO state the order parameter performs one- dimensional 
spatial oscillations, which reduce the pair-breaking effect of the external field. 
Other theoretical studies on the FFLO state include the relation to the mixed 
state H 10 1, the relevance of the shape of the Fermi surface [11, E|, heavy- 



fermion superconductivity [P^| , one-dimensional systems [T3J, the vicinity of 
the tricritical point [15[], and a calculation of the lower critical field |16[]; to 
mention only a few. Recently, it has been shown by one of us [ 17| , that for 



superconductors with cylindrical (circular) Fermi surface one finds, at low 
temperature T, several different types of two-dimensional periodic structures 
with lower free energy than the traditional FFLO-state. Experimentally, 
while indications of spin paramagnetic effects have been observed, no clear 
evidence for an FFLO-like phase has been obtained so far. The recent dis- 
covery of several classes of compounds with nearly isolated conducting layers 
has, however, opened new possibilities to detect the FFLO state. 

In most treatments, either the "pure" orbital effect or the "pure" spin 
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effect has been investigated. Here, we study, following Bulaevskii [18| , a sit- 
uation where both effects must be considered. The relative importance of the 
spin effect is enhanced by means of the following geometrical arrangement: 
Consider a quasi-two-dimensional thin film with a large component Hu of 
the applied magnetic field parallel to the film, and a small perpendicular 
component H±. The smallness of the perpendicular component makes the 
spin effect comparable with the (much stronger) orbital effect. The upper 



critical field for such a situation has first been calculated by Bulaevskii [TJ 
and later studied in the vicinity of the tricritical point | 15| . Still later, Bu- 



laevskii's results were rediscovered in a different theoretical framework, and 
generalized to arbitrary temperatures and d-wave superconductivity by |19 



The essential result is a non-monotonic critical field curve B^iT) which con- 
sists of several pieces belonging to different values of the Landau quantum 
number n. Such a behavior is in contrast to standard Ginzburg-Landau (GL) 
as well as microscopic theory, where the pair-wave function is always found 
to be in the lowest Landau level n = at H c2 . Obviously, the gap at the 
Fermi surface opened by H\\ forces the single electrons in pairing states with 
higher Landau levels n > 0. We note that formally similar upper critical field 
curves may also occur in a different physical context, at extremely high fields 
when Landau quantization of the single electron levels must be taken into ac- 
count PO, ETL E2L E31. For the present situation, Landau quantization of the 



single electron levels can be neglected and the quasiclassical approximation 
holds. 

The present communication addresses the question of the actual structure 
of the spatially inhomogeneous state below some part of the upper critical 
field curve characterized by a Landau quantum number n. This state is of 
interest because it is a consequence of a competition between two different 
pair breaking mechanisms. It is also of practical interest because the portion 
of the H, T plane where it occurs may be easier accessible experimentally than 
the one corresponding to the "pure" FFLO state. Note also, that a perfect 
orientation of the external field is difficult to achieve; as a consequence, a 
small perpendicular component is always present even if one tries to avoid 
it. 



2 Calculation and theoretical background 



We use as a framework the quasiclassical Eilenberger equations [24|, prop- 



erly generalized to include the coupling between electron spins and magnetic 
field |p5 |. The notation used in the present work is the same as in Ref. [26 | . 



We restrict ourselves to pure superconductors with an isotropic gap and a 
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circular Fermi-surface. The quasiclassical equations consist of the transport 
equations for the Greens functions f,f + ,g and the self-consistency equations 
for the order parameter if) and the vector potential A. The transport equa- 
tions are given by: 



2uj s + M F {k) (V R - i(2e/hc)A)] f(R, k, u s ) = 2iJj(R)g(R, k, u 8 ), (1) 
2to s - M F (k) (V fl + z{2e/hc)A)] f + (R, k, cu s ) = 2ip*{R)g(R, k, u a ), (2) 

g(R, k, u s ) = (l - f(R, k, u s )f + (R, k, ^ . (3) 

Here, the Fermi velocity is denoted by v F (k), the Zeeman term is contained 
in the quantity uo s = uj\ + i[iB, where iO\ are Matsubara frequencies, \x is the 
magnetic moment of the electron and B is the magnitude of the induction. 
The self-consistency equations are given by 



N D i 

2nk B Tj2- + log (T/T c 



nk B TY^ 

1=0 



d 2 k' 



V R x (p{R) + AnM(R) 
d 2 k' 



f(R,k',Lu s ) + f(R,k>, 
(4) 



16n 2 ek B TN(E F ) g 



1=0 



Air 



v F (k')Qg(R, k , u> s 



M(R) = 2v 2 N(E F )B(R)-4nk B TN(E F )vJ2 / 



(5) 



(6) 



The magnetization M comes from the spins; the second term in Eq. (|6]) is due 
to the interaction between spins and supercurrent and may be neglected in 
the high /t-limit. These equations have to be supplemented by an expression 
for the Gibbs free energy G. A functional G leading to the above Euler- 
Lagrange equations may be constructed as a straightforward generalization 
of Eilenbergers original expression jP]] and will not be written down here. 
As a first step, the differential operator in the transport equations has 



to be inverted using standard methods p7| . This calculation reproduces the 
equation for B c2 obtained previously |18, [H| , as well as the eigenfunctions of 
the linearized gap equation, which are given by 



<f>n{k, r) 



A 



K± 

-2h iy 



Vk) 



(7) 
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Here, x, y are the coordinates in the superconducting layer, n — 0, 1, 2, . . . 
is a discrete (Landau level) quantum number, k is a continuous quantum 
number, and H n is the Hermite polynomial of order n. The remaining quan- 
tities in (0) are given by k± = — £?sin0 and = fc//cj_, where 9 is the 
(small) angle between the external magnetic field and the superconducting 
plane. The value of A may be chosen to normalize the complete, orthogonal 
set of eigenfunctions (0). The same form of the normal state vector potential 



as in previous work |T9| , |28| has been chosen. 

We want to study the structure of the (stable) inhomogeneous state below 
B c2 for arbitrary n. For given temperature T, angle G, and ratio T c /Ep, the 
correct value of n is obtained by solving numerically the equation for B C 2 |19| . 
Given this n, the order parameter below the corresponding branch of the B C 2 
curve can be constructed as a linear combination of eigenfunctions belonging 
to the infinitely degenerate set (0), which are labeled by the continuous 
quantum number k. We construct a general state with an order parameter 
which is quasi-periodic (periodic up to a phase factor) in two dimensions. 
Abrikosov's method [p{|, used originally for n = and T ~ T c , may be 
applied in a straightforward manner to the present case of arbitrary n. It 
leads to an order parameter ip n which is given by 



a 2 

h 



■ exp ( —k m {X + Y cos a) ) h (sin a (Y — mb)) , (8) 



where 




Here, a is the angle between the two primitive lattice vectors a, h (of length 
a, b) spanning the unit cell of the periodic structure we are interested in. 
An oblique system of coordinates X, Y has been used in (§) whose relation 
to Cartesian coordinates x, y is given by x = X + Ycosa, y = Ysina. In 
the course of the calculation leading to the continuous numbers k have 
been restricted to a discrete set k m = {2%h/a)m with integer m. The flux 
due to the perpendicular component B± = B sin 9 must be quantized; the 
corresponding fluxoid quantization condition has been used in the derivation 
of (|8|). It takes the form n±ab sin a = 2irh if we assume that each unit cell 
carries one flux quantum. Periodic order parameters for arbitrary n have 
been reported previously |5D], [HJ] in the literature. 

Among all possible structures @ the stable one, i.e. the one with the 
lowest free energy, must be found. To perform this task the quasiclassical 
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equations have to be solved near B C 2 and the free energy has to be calculated 
up to terms of fourth order in the order parameter magnitude. In contrast 
to the GL region, the coherence lenght is finite at low T; therefore all pow- 
ers of the order parameter gradient must be taken into account. This large 
calculation, which generalizes Abrikosov's solution to arbitrary T, has first 
been performed by Eilenberger [32] for the ordinary vortex lattice, i.e. for 
n = 0. Later, Eilenbergers mixed state calculations have been generalized by 
Rammer and Pesch |33| to include strong coupling effects. The present calcu- 
lation generalizes Eilenbergers work to arbitrary Landau quantum numbers 
n. It makes extensive use of calculational techniques ]33], ^7|, [32], |35| de- 
veloped by previous workers in this field and comprises the following main 
steps: (1) Solution of the linearized transport equations, (2) expansion of 
transport equations and free energy for small order parameter, (3) solving 
Maxwells equation, and (4) performing a large number of spatial and mo- 
mentum integrations. While all second order free energy contributions can 
be calculated exactly, one has to perform an asymptotic approximation in 
integrals appearing in fourth order terms in order to obtain a simple final 
result. This approximation is similar to the one suggested by Delrieu |35[ 
and holds for not too low T (see the discussion following Eq.(6) of Delrieu's 
paper fl35| ). 

After a long calculation, which will be described in more detail elsewhere, 
one obtains the following expansion (note that dimensionless units 
be used in what follows) for the free energy G: 



281 will 



G = G 



a 



(10) 



Here, a denotes the spatial average of tp , our small expansion parameter. 
The term G stems from the spatially constant part of the magnetic induc- 
tion [32]; it is unimportant in the present context since it does not depend 
on the structure variables a, a, b. The second order contribution is given by 



r°° i 

lnt + t / ds — — 
Jo smh( 



h(st) 



l-cos^Rs)/^ 2 ^) 



'1L 



where t = T/T c , and B, B± are spatially averaged values of B(R) and B±(R) 
respectively. The function fi depends on the Landau quantum number n and 
is given by 



= exp( 



x , 



LJx) 



(12) 



where L n is the Laguerre polynomial of order n. The second order term ([TT]) 
is also independent of a, a, b, but the equation G^ 2 ^ = has to be solved in 
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order to find the upper critical field and the corresponding value of n (this 
equation for B c2 agrees with the one given in [|l9|j ). 

The new information is contained in the fourth order term The 
magnetic part, say G$ , of is of order 1/k 2 . This magnetic contribu- 
tion collects all the terms related to Maxwells equation. If Eilenbergers, 
or Abrikosovs results are to be rederived, as a special case of the present 
theory for H\\ = {n — 0), all contributions to the free energy, including 
G^, have, of course, to be taken into account (one can successfully check 
Eq.([L(]) against these classical results; the same is true for the FFLO upper 
critical field which may be derived by solving G^ = in the limit Bj_ = 
). In the present communication we consider the high-K limit and neglect 
the contribution G^/ . This approximation is used by most workers on the 
FFLO state, in particular in Larkin and Ovchinnikov's original paper ||. 
A complete treatment including this magnetic contribution as well as mag- 
netization curves for the present geometry, will be reported in a separate 
publication |5B[. Then, the final result for the fourth order term is given by 

-X +oo +oo 

g {4) = - a E E fiMSu, (is) 

l=— co j=— oo 

where fi has been defined in Eq.(^). The Matsubara sum Sij is defined by 

bl <i ~ 1^ — ; : 7s72 ' 

i=o oof [oof + B ±XlJ /2) 

and Xij is given by 

7r 2 /L\ 2 l9 /a\ 2 9 , cos a . , 

*W = T-(-) 1 + T f-27rlj — . 15 
sm a \a/ \LJ sma 

The minimization with respect to the flux-line structure is done at fixed 
B±, B\\ (the minimization with respect to B±, B\\, a has to be done in separate 
steps). Therefore, stability conditions should not be applied with respect to 
the variables a/b, a but with respect to a/L, a, where L is a magnetic length 
defined by L 2 = 2/B±. Then, using the results for a/L, a, the ratio a/b is 
finally determined with the help of the fluxoid quantization condition. 

The fourth order term ( |T3| ) has a simple interpretation. As is well known, 
near T c the conventional flux-line structure is completely determined by a 
single quantity, Abrikosov's geometrical factor Pa = ip 4 /ip 2 (the bar de- 
notes a spatial average). If this quantity is computed using a general order 
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parameter ip n , as given by Eq. (|8|), one obtains, after a lengthy calculation 

+00 +00 

Pa,u= E E ftM- (is) 

!=— 00 j=—oo 

One notices that the free energy ( |I3"D differs from /3^„ only by the temperature- 
dependent factor Sij modifying each term in the double sum (for n — 0, fi^ n 
takes its lowest value 1.1596 for a triangular lattice). The factor B±_xij/2 in 
Si j, (see Eq. (|14|), vanishes in the GL-limit; it reflects nonlocal correlations 
in Eilenbergers transport equations and implies an additional "microscopic" 
temperature dependence of the stable magnetic flux structure. One expects, 
that for higher T, the GL approximation (3a may be used as a first approx- 
imation for the complete free energy (|13|). It should also be mentioned that 
the electrons magnetic moment /i, which is the typical spin pair breaking 
parameter, does not explicitly occur in Eq. flTJ) (but only in G, G (2) , and 
Gm). The spin effect influences, however, the free energy - and the resulting 
magnetic flux structure - in a very decisive way since it determines the value 
of the Landau quantum number n in Eq. ([13D . 



3 Results and discussion 

To calculate numerically the stable flux structure the following two steps 
have to be performed: (1) The equation G^ = has to be solved for a 
given temperature T and angle O to obtain the upper critical field and the 
corresponding value of n. (2) Using these values the stable flux structure 
has to be determined by calculating the absolute minimum of ([U|) with 
respect to a/L, a. Throughout the numerical calculations a value of fi = 0.1 
(/i = (n/2)kBT c /Ef in the present system of units) has been used. The 
calculations have been performed (for various n in the range between and 
30) at two different temperatures t = T/T c = 0.5 and t = 0.2, in order to 
investigate the influence of the above discussed low-temperature correlations. 



The minimization procedure turned out to be more involved |37j than one 
would expect from the simple form of Eq. (|T3|). One finds a large number 
of minima belonging either to different lattice structures or to different basis 
vectors of the same lattice. The number of minima increases rapidly with 
increasing Landau quantum number n. The quantity shown in all plots is 
actually G$ , which is G^ multiplied by a factor of 16t 3 /7x(3) ~ 2t 3 (in 
order to normalize it for t — > 1 to the GL result). Fig. 1 shows a contour 
plot for n = 10, t — 0.5 of G^ in the most important part of the a,a/L- 
plane; higher values of a/L need not be considered, since they correspond 
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to lattices with exchanged axes a, b. The minima of the free energy are 
visible as dark regions. Two different types of minima may be distinguished 
in Fig.l (and similar plots). A minimum of the first type is represented by 
an isolated dark region and corresponds to an "ordinary" two-dimensional 
lattice (several equivalent points in the a,a/L-plane, which have the same 
free energy and belong to the same lattice, but to different basis vectors, 
may be found in Fig.l). The second type of minimum can be found in the 
straight (a-independent) dark regions in the lower part of Fig 1. Near these 
"flat valleys" the free energy depends mainly on a/L, a small dependence on 
a along some lines of fixed a/L exists, but is invisible in Fig 1. The minima 
in these flat valleys will be referred to as "quasi-one-dimensional" lattices, 
because of the small energy barrier along the lines of constant a/L; a more 
detailed discussion will be given below. 

The phases with n — 0, 1, 2 occupy a relatively large portion of the phase 
diagram [18|, [19| and are therefore of particular interest in view of a possible 
experimental detection of these new states. For the lowest Landau index 
n = 0, we recovered at both temperatures, as expected, the conventional 
triangular vortex lattice as stable state (global minimum of the free energy). 
To create a state with n = 1, a possible choice for the tilting angle of the 
external magnetic field is = 0.3 (degree) at t = 0.5; another possible 
choice is 6 = 1.2 at t = 0.2. The stable structure found for n = 1 is, at 
both temperatures, of the quasi-one-dimensional type. A contour plot of the 
magnitude of the order parameter is shown in Fig. 2; the structure consists 
of chains of vortices separated by a single line of vanishing order parameter. 
The latter may be considered as a fragment of the one-dimensional FFLO- 
state; for higher n local minima appear whose order parameters look similar 
to Fig. 2 except that n lines of order parameter zeros appear between the 
vortex chains. 

Let us discuss these quasi-one-dimensional structures - which have been 
verified as local, but not necessarily global, minima for many higher Landau 
indices in the range n < 30 - in more detail. They consist of two subsys- 
tems, vortices and one-dimensional FFLO-type oscillations. Both of these 
are of the "ordinary" type, i.e. they are basically unaffected by their mutual 
interaction; the phase change when encircling a vortex is 2tt and the wave 
length as well as the behavior of the phase of the one-dimensional periodic 
structure agree with that of the one-dimensional FFLO state. For all of these 
states, the coupling energy between neighboring vortex chains is small. As 
a consequence, the energy barrier preventing motion of the vortex chains in 
the direction parallel to the chains is very small (for n = 1 the a-dependence 
shows up in the sixth digit of the free energy). Thus, even small fluctuations 
will destroy the periodicity in the direction parallel to the chains, making the 



9 



material effectively one-dimensional. These decoupling of the vortex chains 
leads to a quasi- continuum of unit-cells (similar to the one shown in the lower 
part of Fig. 1 for n = 10) which is located near a/L = 1.03 for n — 1. Ex- 
perimentally, one expects a pronounced change in transport properties if, by 
decreasing the angle 9, the system makes a transition from n — to n — 1. 

As regards the next higher Landau level index, n = 2, stability calcula- 
tions have been performed for the two sets of parameters = 0.17, t = 0.5 
and O = 0.7, t = 0.2. At t = 0.5 a triangular lattice is realized as the 
state of lowest free energy. At t = 0.2 a slightly distorted triangular lattice 
(a = 60.05) is realized. The distortions are due to the microscopic factor Sij 
in Eq. (0). Thus, the influence of the microscopic correlations on these low- 
n states is visible but very small. The order parameter for n = 2 is shown in 
Fig. 3. The triangular unit cell contains three order parameter zeros. These 
zeros have different vorticity; the phase change when encircling a zero is +27T 
for two of the vortices and — 27r for the third. Recall that, for arbitrary n, 
the unit cell as a whole carries exactly one flux quantum. 

For n = 3 at t = 0.5 we found again a quasi-one- dimensional stable 
state. The order parameter for n = 3 looks similar to Fig. 2, except for 
three dark lines appearing between the vortex chains. No stable quasi-one- 
dimensional states have been found for higher n (n > 3). In the interval 
3 < n < 10 the global free energy minima are given by various types of two- 
dimensional lattices, including triangular, quadratic and oblique structures. 
With increasing n the influence of the microscopic correlations on the lattice 
structure increases; the stable state does not only depend on n but may 
change between t — 0.5 and t — 0.2 (keeping n constant). 

The phases with higher n correspond to very small values of B (e.g. for 
n = 8 and t = 0.5 one has = 0.05). One expects the resulting stable states 
to be already very similar to the FFLO state, which corresponds to the limit 
n — > oo, = (B± = 0). This FFLO-limit n — > oo of the present theory is a 
subtle step. An analytical limit of the calculation cannot be performed since 
in the present formulation extensive use has been made of the fluxoid quanti- 
zation condition, which becomes meaningless if no perpendicular component 
exists. Nevertheless, one expects for general reasons of continuity, that the 
present theory comprises the FFLO theory in the sense, that, with increasing 
n, the numerical results become more and more similar to the FFLO state. 
Thus, for t = 0.5 we expected, as a natural candidate for a smooth transition 
to the FFLO state, the quasi-one-dimensional state to be stable for large n. 
At t = 0.2 we expected, in view of JL7|, a two-dimensional structure to be 
stable for large n. 

However, contrary to our expectation, no single stable lattice could be 
found in our calculations which have been performed up to n = 30 for t = 0.5 
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and n = 20 for t = 0.2. Contour plots of the free energy look similar 
to Fig. 1 but with an increased number of local minima at higher n. At 
both temperatures the transition from n to n+1 (keeping t fixed and varying 
G) leads generally not to the same but to a different lattice; among these 
stable states we found all kinds of oblique, triangular and rectangular two- 
dimensional structures. Quasi-one dimensional states, which are stable at 
n=l,3, have also been found as local minima in these high-n calculations, but 
had always higher energy than the minima corresponding to two-dimensional 
states. The difference in free energy of the different local minima decreases 
with increasing n. 

At present, we are unable to decide whether this strange vortex-liquid 
like limiting behavior, which would imply the instability of the FFLO state, 
is an artifact of our calculations or corresponds to reality. There are some 
arguments in favor of the latter possibility: Adding a perpendicular compo- 
nent of the magnetic field, no matter how small, implies flux quantization; 
this changes the conditions which determine the stable state in a discon- 
tinuous way. Further, the size of the unit cell and the number of vortices 
(with different vorticity) per unit cell increases with n. It seems reasonable 
to assume that the energy differences between different vortex arrangements 
(lattice structures) become small for large n. More work is required to settle 
this question. 

4 Conclusion 

We considered a geometrical arrangement with competition of orbital and 
spin pair breaking effects in a two-dimensional superconducting state. This 
arrangement comprises the traditional vortex state (for B\\ = 0) as well as 
the FFLO state (for B±_ = 0). Minimizing the free energy we found several 
new structures below B c2 with pairing states in Landau levels n > 0. These 
include two-dimensional structures with vortices of different vorticity as well 
as mixtures of one-dimensional periodic order parameter oscillations and vor- 
tex chains. The latter show a fluid-like behavior, as regards vortex motion in 
the direction parallel to the chains. This feature could be used, besides more 
local spectral techniques, to identify this new state experimentally. Topics to 
be treated in future work include the high-n limit, contributions from finite 
values of k, magnetization curves, a more systematic account of unit cell 
structures, and the transitions between pairing states belonging to different 
n. 
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FIGURE CAPTIONS 



1. Contour plot of the fourth order free energy contribution G N for n = 
10, t = 0.5 as a function of a,a/L. Minima of free energy must be 
sought in dark regions. 

2. Contour plot of as a function of x,y. Stable structure for Landau 
level n — 1 at t — 0.5 with lattice parameters a/b = 2.75, a = 69.83. 
A small energy barrier exists against translation of the vortex chains 
relative to each other. 

3. Contour plot of \ip\ 2 as a function of x,y. Stable structure for Landau 
level n = 2 at t = 0.5 with lattice parameters a/b = 1.00, a = 60.00 
(a/L = 1.90). The unit cell carries one flux quantum and contains 
three order parameter zeros. 
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Fig. 3 Klein,Raiiier,Shimaliara 



